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ABSTRACT 

The radiative and jet efficiencies of thin magnetized accretion disks around black holes 
(BHs) are affected by BH spin and the presence of a magnetic field that, when strong, could 
lead to large deviations from Novikov-Thorne (NT) thin disk theory. To seek the maximum 
deviations, we perform general relativistic magnetohydrodynamic (GRMHD) simulations of 
radiatively efficient thin (half-height H to radius R of HjR ^ 0.10) disks around moderately 
rotating BHs with a/M = 0.5. First, our simulations, each evolved for more than 70,000rg/c 
(gravitational radius Pg and speed of light c), show that large-scale magnetic field readily 
accretes inward even through our thin disk and builds-up to the magnetically-arrested disk 
(MAD) state. Second, our simulations of thin MADs show the disk achieves a radiative ef¬ 
ficiency of T/r ai 15% (after estimating photon capture), which is about twice the NT value 
of ? 7 r ~ 8% for a/M = 0.5 and gives the same luminosity as a NT disk with ajM 0.9. 
Compared to prior simulations with < 10% deviations, our result of an ^ 80% deviation sets 
a new benchmark. Building on prior work, we are now able to complete an important scaling 
law which suggest that observed jet quenching in the high-soft state in BH X-ray binaries is 
consistent with an ever-present MAD state with a weak yet sustained jet. 

Key words: accretion, accretion discs, black hole physics, hydrodynamics, (magnetohydro¬ 
dynamics) MHD, methods; numerical, gravitation 


1 INTRODUCTION 

Black hole (BH) accretion systems can operate as efficient engines 
converting gravitational potential energy and BH spin energy into 
radiation. Quasars, active galactic nuclei (AGN), X-ray binaries, 
gamma-ray bursts, and other BH accretion systems enter a thin ra¬ 
diatively efficient mode when the luminosity is between ~ 10“^ and 
~ 0.5 times the Eddington luminosity l Abramowicz et ar] |1995[ 
|Narayan & Yi|1995) . The |Novikov & Thome P973 i, NT, calcula¬ 
tion for the thin disk radiative efficiency assumes emission termi¬ 
nates at the inner-most stable circular orbit (ISCO). However, ac¬ 
cretion is driven by magnetic stresses via the magneto-rotational 
instability (MRI) ( [Balbus & Hawleyj|1991| |1998a|l or magnetic 
braking by large scale field threading the disk (Blandford & Payr^ 
|1982| > or BH jBlandford & Znajek|1977[ l, hereafter BZ. A funda¬ 
mental question in accretion theory has been whether the NT model 
is applicable to actual magnetized disks where magnetic flux cou¬ 
ples to a rotating BH and increases the radiative efficiency near the 
ISCO ( |Gammie|1999[[&r3lik|1999^ . Any extra radiative efficiency 
could, in principle, affect BH spin estimates obtained from Soltan’s 
argument (e.g., see [Shapiro |20()^ or from co ntinuum emission 
spectrum observations jMcClintock et al.|201 ijl. 


An upper limit to the amount of magnetic flux a disk and BH in 
steady-state can support is achieved when an accretion flow reaches 
the so-called magnetically-arrested disk (MAD) state, when mag¬ 
netic flux accretes until magnetic forces pushing out balance gas 
forces pushing in (jlgumenshchev et al.|2003||Narayan et al.|2003[ 


Igumenshchev|2008[ [Tchekhovskoy et al.||2011 [[McKinney et al.| 


2012 1. If Jet power scales as Pj oc for magnetic field B and 


BH spin a/M ( [Blandford & Znaj^|1977| >, then the MAD state 
leads to the maximum jet power. The existence of such a maxi¬ 
mum well-defined magnetic flux for a given spin is assumed when 
using transient jets in BH X-ray binaries to test the BZ prediction 
of power vs. spin ( jMcClintock et al.|2014^ . Destruction of such a 
strong magnetic flux could lead to such transient jets without the 
need for BH spin Pgumenshchev|200^ [Dexter et al.|2014| (, while 
BH spin measurements may be more reliable for thin MADs that 
could have angular momenta aligned with the BH spin axis near the 
BH ( [Polko & McKinney[2015[ (. 

General relativistic magnetohydrodynamic (GRMHD) simu¬ 
lations of thin disks have found different degrees of deviations from 
NT theory, but in these cases there are no more than 10% fractional 


deviations for the radiative efficiency (Shafee et al. 

20081 [Noble[ 

[et al.|2009[|Penna et al.|20101|Noble et al.|2010||2011 

. Differences 
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in results are likely due to the amount of magnetic flux that threads 
the BH and disk ([Penna et al.[2010[(. Unlike other choices for initial 
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conditions, the amount of magnetic flux threading the BH and disk 
in a MAD state is independent of the initial magnetic field strength 
or geometry as long as there is a plentiful supply of magnetic flux 
( [McKinney et al.|20T2{|Tchekhovskoy & McKinney|2012^ . Numer¬ 
ous non-MAD GRMHD simulations of relatively thick disks have 
been performed Wliers&Hawley|2003JJMcKinney&^^ 


miej2004|^ De^^llierertay200^ Mcffinneyj2005 J200^ Komis 


sarov & McKinney|2007[|McKinney & Narayan|2007[ Fragile et al. 
2007[|McKinney & Blandford|2009[|Narayan et al.|2012^ . Simula- 
tions of MADs have so-far involved only relatively thick disks with 
no cooling ( jTchekhovskoy et al, [20 11 [[McKinney et al.\20i2) and 
thick disks with radiative transfer ( [McKinney et al.[2015| l. 

We perform fully 3D GRMHD simulations of radiatively ef¬ 
ficient thin (half-height H to radius R of H/R x 0.1) disks around 
rotating BHs (dimensionless spin, ajM = 0.5) that reach the MAD 
state. We expect the MAD to maximize the radiative efficiency of 
thin disks and lead to the maximum deviations from the NT effi¬ 
ciency. First, we test how easily magnetic flux advects inward by 
threading the disk with weak magnetic field. Second, we setup a 
disk that is initially nearly MAD, let it become MAD to large ra¬ 
dius, and then measure the radiative efficiency. 

We describe the physical and numerical setup in ® present 
results and provide discussions in ^ and summarize in 


2 FULLY 3D GRMHD THIN MAD MODEL 

For this study we use the HARM GRMHD code ( [Gammie et al.[ 
2003|l with variou s improvements ( [McKinney & Blandford|2009 
McKinney et al.[f2012^ with two physical setups, both with spin 
a/M = 0.5. Each was run with a low resolution, medium resolution, 
and high resolution for convergence testing and the high resolution 
was used for detailed analysis. We name the simulations in this 
work MADxy where x=i,f signifies whether the disk starts with 
enough flux to be nearly MAD in the initial, i, or only MAD in the 
final, f, state. y=HR,MR distinguishes the grid as high or medium 
resolution. 


of magnetic flux threading the BH and disk, we provide a large sup¬ 
ply, and the amount of supply ends up not mattering because it is 
vast compared to what ends up on the BH or in the disk near the 
BH. 

The coordinate basis (/i-component of the vector potential is 
matched to the radial gas pressure profile so that the plasma param¬ 
eter p = Pgisll'b (Pb as magnetic pressure), and thus the number of 
critical MRI wavelengths spanning a full scale-height of the disk, 
is a tunable constant for all radii. For MADiy runs, we found it use¬ 
ful to transition from MAD-flux levels to a weaker sub-MAD field 
beyond a transition radius, tq.mad = 30, with a transition function 
Fj, = 0.5 ± ^ arctan((r - ro,MAD)/4). The vector potential is then 
calculated by radially integrating 


{sF- + F+) sin 9'' ■ 


K 


(- 0 . 6 / 2 + 3 / 2 )/,^ 


if r > r„ 
if r <= Ttr 


( 1 ) 


where is the field threading the disk mid-plane at time t=0 nec¬ 
essary to start with s = 2 and p = 200/= 50 inside ro MADi and 
with p = 200 elsewhere. Thus, the disk within tq^mad initially has 
5d ~ 0.7, where 5d is the number of MRI wavelengths within two 
disk scale heights, Sd £ 0.5 where the MRI is suppressed jBalbus 


& Hawley 


mney et a 


1998b , and MADs typically have had Sd £ 0.25 ( [McK- 
.[2012 I. This means the MRI is marginally suppressed 


in MADs. The field threading the BH horizon and within the ISCO 
matches to the disk values at a transition radius of rtr = 12. /r = 4 
ensures the initial flux residing within is small enough (« 2.5 
times smaller than what ends up on the BH plus out to r,r) that the 
BH and the region between the horizon and r„ can slowly build-up 
to a MAD level. 

To ensure a MAD can readily develop for a thin disk, our 
MADf model starts with sub-MAD flux at t=0, with p = 200 out¬ 
side r,r. The solution for is the same as the r > r,j branch of 
Eq.[T]but with ^o mad and a linear taper in the field from r,, 

to the ISCO so there is no initial magnetic flux any closer. 


2.1 Physical Model 

The initial disk is Keplerian with rest-mass density that is Gaus¬ 
sian in angle with a half-height-to-radius ratio FljR 0.1 and a 
power-law radial distribution p oc r^-^. Since this is not an equi¬ 
librium solution near and inside the ISCO at rjsco. the density is 
tapered as p —> p(l - .95(risco/R)^^) and truncated at nsco- The 
gas internal energy density Cgas is obtained from the vertical hy¬ 
drostatic equilibrium condition of HjR x h = c^/Ujo, for sound 
speed Cj s; -\jrPg.^slp and rotational speed Vmt chosen to be Keple¬ 
rian with velocity Vjc = (rlrg)l{(r/rg)^F + ajM). In order to seed 
the MRI we use a 10% random perturbation to the pressure that 
obeys the ideal gas law Bgas = (F - IfCgas- We choose F = 4/3 
because the disk is expected to be radiation dominated near the 
BH. We set an atmosphere for the disk with p = 10“'*(r/rj)“^^^ and 

Cgas = 10-^(r/rg)-5/2 

We choose an initially purely poloidal large-scale magnetic 
field. With other low resolution models already performed, we 
know ahead of time the magnetic flux reached on the BH and 
throughout the disk in the quasi-steady MAD state, so we then tune 
our initial conditions to be either nearly MAD (for MADi) or very 
sub-MAD (for MADf). Net magnetic flux on the BH and through 
the disk out to some radius is conserved and only changes by accre¬ 
tion or expulsion through that radius. To reach the maximal amount 


2.2 Numerical Grid and Density Floors 

We use similar numerical grid mapping and identical boundary 
conditions as [McKinney et al!]j2012[ l, except a smooth arctan tran¬ 
sition from exponential to hyper-exponential radial grid spacing 
is used. For our high (mid) resolution grid we choose resolution 
N,. X Ng X = 192x96x208 (168x64x152) with radial grid span¬ 
ning Bin ~ OJSrn (horizon radius th) to Bout = 

The 0-grid spans from 0 to tt with the mapping of [McKinne y[ 
[et al.[ ( [2012^ , but with ujet = 0 used for our entire disk since it is 
thin at all radii. We adjust the constants of the mapping so that at 
rjrg = 6,30,100 we resolve the disk vertically with 56,62,60 cells 
within ±2H= 2 • 0. IR. We flare the inner radial grid to avoid unnec¬ 
essarily small limiting time-steps. The polar boundary conditions 
are transmissive (see appendix in [McKinney et al.|2012[ l, and the 
(/i-grid is equally spaced from 0 to 2n with periodic boundary con¬ 
ditions. 

Apart from our direct resolution studies, we also measure con¬ 
vergence quality factors for the MRI and turbulence correlations 
lengths. Our values demonstrate good resolution and are reported 
in Section iTS] 
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Efficiency of Thin MADs 3 



Figure 1. Snapshots of the initially sub-MAD (MADfHR - left panel) and nearly-MAD (MADiHR - right panel) models at t=0, showing p in color (with 
legend) with magnetic field lines (black lines). In MADfHR there is no magnetic field initially present within the ISCO or on the horizon and the field is weak 
througout the disk. On the other hand, in MADiHR, a sufficient magnetic flux is present throughout the disk and on the horizon so that the disk will quickly 
become MAD out to large radius. The higher density surrounding the horizon in MADiHR at t=0 is the result of the density floor activated due to the presence 
of an initially strong magnetic field. 


2.3 Cooling 


To keep the disk thin and radiatively efficient, we implement an ad 
hoc cooling function similar to|Noble et al.|pOT^ with 


Sv = 



( 2 ) 


as a source term to the conservation equations, = 5v, for gas 
stress-energy tensor T and 4-velocity u,. This assumes an isotropic 
comoving rate of thermal energy loss dUldr. dU is computed so 
that cooling decreases internal energy toward a target temperature 
set by the desired H/R x: 0.1 and local density and radius. 

Radiative cooling is efficient for a NT thin disk with a cooling 
timescale on the order of the orbital timescale jNovikov & Thoni^ 
|1973| >. The NT model assumes radiation gets out only by diffu¬ 
sion, while our disk is clumpy and can be optically thin near the 
ISCO ( |Zhu et al.|2012) . Radiative GRMHD MAD simulations are 
required for an accurate cooling timescale, which could be faster 
than for literal NT disks ( [McKinney et al.|[2015^ . We only try to 
keep the disk thin, and we settled upon a graded cooling timescale 
dr so that 


dU 

dr 


pCtarget 


dr 


— pTtarget 


^^k/'^cooI 

(i^) ^ 

VA\f 

H 

Tcooi/Hk > > 2A? 

1 

. 2Ar 

(4L)<2Af 

V VAlf / 


(3) 

where VAif is the Alfvenic frequency and thus HIvau is the Alfven 
time across one disk scale height, 2At is the code time-step with a 
prefactor 2 for stability, and Dk = The target temperature is 
from thin disk theory and is isotropic across cylindrical radii R, 


0.1 R 


where the denominator uses spherical radius r to control the tem¬ 
perature along the polar axes. If Tgas = Pgas/p < Ttarget then 
dU/dT = 0. To ensure there is always rapid enough cooling, we set 
Tcooi = 0.1 • 2;r that enforces cooling on 1/10* the orbital timescale. 


We also use a temperature ceiling that acts as another cool¬ 
ing term, operating after each timestep. We limit the temperature 
where /? < 0.1 and Tgas > Ttarget- Unfortunately these losses were 
not tracked in the version of HARM used. When we sum the ra¬ 
diative losses from the cooling function during our analysis, in or¬ 
der to compute the contribution lost to the temperature ceiling, we 
restarted MADiHR at a late time with the ceiling turned off. The 
additional heat is then captured by the cooling function which we 
can measure for the radii near the BH where the temperature ceil¬ 
ing mattered. Accounting for all cooling effects gives us a nearly 
constant radial profile for total efficiency, and so the procedure cor¬ 
rectly recovers how a steady-state flow should behave. 

Rest-mass density is injected whenever TP^/p > 200 in order 
to ensure the code remains stable, and internal energy density is 
added when figas/p > 100 (though this limit is never reached). 

As seen in Fig. the cooling keeps the disk very close to the 
target scale height. After a very short (< 5000rg/c) initial period of 
disk evolution, the disk height is very stable. 


2.4 Diagnostics 

The disk’s geometric half-angular thickness (H) per radius (R) is 

{lp(e-dordAg^) 


\l/n 


where 




do(r,(l>) = 2 * 


(fpdAg^) 

(fgp(S - ^l2)"dAe^'j 
i/iT 


l/n 


(4) 


(5) 


and with surface differential dAg^. Typically we use n = 1 but also 
compare to using other n. 

The radiative energy and angular momentum lost are 

E(r) = r KiAdAg^ = - -a: r r (6) 

J JDt JRc&p 
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Figure 2. Evolved snapshot of the MADfHR model at f « 70, GOOr^/c showing log of rest-mass density in color (with legend) in both the z - x plane at 
y = 0 (top-left panel) and y - x plane at z = 0 (top-right panel). Black lines with arrows trace field lines, where thicker black lines show where field is 
lightly mass-loaded. The bottom panel has 3 subpanels. The top subpanel shows the mass accretion rate through the BH (Mr)- The middle subpanel shows 
the magnetic flux threading the BH (Tr), normalized so that order unity is a dynamically-substantial amount of magnetic flux. The bottom subpanel shows 
the gas efficiency measured on the horizon (? 7 R,g) and jet efficiency measured at r = 50r^ (rjj). In this model, the disk initially has a weak large-scale poloidal 
magnetic field, which is transported by the MRI and other magnetic stresses toward the BH. Eventually, by r ~ 55000rg/c, the magnetic flux saturates on the 
BH and cannot grow despite plentiful magnetic flux outside. The MAD state that develops consists of a highly non-axisymmetric clumpy thin disk with a weak 
(’ll ~ l%)jet. 


respectively, for v = for some period Dt over all angles with 

photon capture radius (see The radiation contribution from 

beyond the inflow equilibrium radius (~ for our MADiHR 
model) is conservatively estimated by using NT’s radial depen¬ 
dence. 

The mass accretion rate, energy efficiency, and specific angu¬ 
lar momentum accreted are, respectively. 


M 

'7 

J 



f(r[+pu'-+R^)dAg^ 


[V]// 


(7) 

( 8 ) 

(9) 


where R is the radiation stress-energy tensor and [M}h is the time- 
averaged M on the horizon. These quantities are calculated as totals 
and for the jet and wind. The jet is defined as where 2Pi,lp > 1, and 
the wind is where IPtjp < 1 with an outgoing flow. The radiative 
efficiency and specific angular momentum vs. radius are, respec¬ 
tively, 77 , = -L,(r)l[M]H and = L^(r)l[M]H. 


xiic; u.iiiic;iisnjiiic/ss iiiagiic-tiL. iiua, iiiai iin 

of the magnetic field threading the BH or disk, 

/ dA,^Q.5W\ 


Y(r) ^ 0.74 




( 10 ) 


for radial magnetic field strength B' in Heaviside-Lorentz units, 
where Yh is T measured on the horizon and 


'T(r) = J dA,^0.5IB'l 


i the radial absolute flux passing through a given radius. 

Viscous theory gives a GR ff-viscosity estimate for Vr of 

rri / n\2 L. I _..i J-inn 


■Ga{HIR)^\v^\, as defined in McKinney et al. (2012i, 


Up vise ^ tiij vikxiiiikxVi. Ill ivxv/iviiiiiV'^ wi cii • y 

leading to a definition for the measure of effective a-yiscosity 

Vr 

( 11 ) 

t^visc/ tt" 

whereas yiscosity from local shear-stress in the disk is 

brb^ 

Q'T, ~ rtmag “ . (12) 

Pb 

where we haye dropped other (e.g. Reynolds) terms from the full a 
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Efficiency of Thin MADs 5 


because they are negligible compared to the magnetic term in the 
MAD regions of our simulations. These a’s (see Fig. are ob¬ 
tained by separately volume-averaging the numerator and denomi¬ 
nator in 9, <p for each r with a weight of p and only including mate¬ 
rial with E!Pa < 1 in order to focus the measurement on the disk 
material and not the lower density corona. The quantity mea¬ 
sures the magnetic stress within the dense regions, while the 
measures the actual effective a based upon the radial velocity of the 
dense regions. For example, very little magnetic flux may thread the 
dense regions, but external magnetic torques may push and drive 
angular momentum transport and inflow of the dense regions. ai,,ef! 
is only accurate far outside nsco because of the G term. 


3 RESULTS AND DISCUSSION 

With our fully global 3D GRMHD simulations we can demonstrate 
whether the non-linear MRI and other MHD physics in thin accre¬ 
tion flows is able to transport large-scale magnetic flux towards the 
BH and build-up the magnetic field to the MAD state. We first con¬ 
sider a disk with a weak magnetic field using our MADiHR setup 
to see if magnetic flux transports inward. Then, having confirmed 
thin disks transport magnetic flux efficiently to the BH and inner- 
radial disk, we use our second setup MADiHR to study a MAD 
where magnetic saturation is reached to large radii. This two-step 
procedure is necessary due to the extreme computational resources 
that would otherwise be necessary to build-up large-scale flux to 
the MAD state out to large radii. 

3.1 Demonstration of flux accumulation to MAD state 

First, we consider a thin radiatively efficient disk threaded by weak 
large-scale magnetic field. This model, MADfHR/MR, started in 
the sub-MAD state with s; 200 (a field strength that is ~ 7 times 
weaker than when the MRI is suppressed and ~ 40 times weaker 
than the final MAD state) so that the MRI operates over long times. 
In MADfHR the MRI leads to accretion of magnetic flux inward, 
such that over t x 70, OGOr^/c the BH and disk become MAD out 
to r ~ 15rj. Fig. shows a steady accumulation of flux during the 
first ~ 30000rj/c. The level of flux then remains constant on the 
horizon for another ~ 20000rg/c due to an excess of rest-mass ac¬ 
cumulating in the inner disk that slowly accretes. This denser torus 
forms as the magnetic flux threading the inner section of the disk 
accretes onto the BH during the early evolution and temporarily 
starves some disk material of strong field. The local viscosity is 
then suppressed until magnetic flux from further out in the disk has 
time to transport inward. Once this initial dense torus accretes, the 
disk enters a MAD regime of accretion. 

In the much longer run, MADfMR becomes MAD out to a 
radius of r ~ ISr^ by the time t x 183, OOOrg/c. Plenty more mag¬ 
netic flux exists at larger radii, but it would presumably contribute 
to the MAD region of the disk at later times. The physics of mag¬ 
netic flux accumulation involves the weak field MRI driving accre¬ 
tion, magnetic braking via the wind transporting material and mag¬ 
netic flux inward, and may involve the field penetrating coronal-like 
material in the disk (e.g., see|Rothstein & Lovelace|2008||Beckwith| 
|et al.|2009) . At early times, magnetic flux at higher latitudes trans¬ 
port more readily than those in the equator, but eventually magnetic 
flux threading the disk is also transported to the region close to the 
BH. Whatever is the physics behind the saturation of magnetic flux, 
it has reached a maximum far below the supply, yet far beyond the 
initial value. This demonstrates that large-scale magnetic flux can 


readily accumulate even with H/R 0.1 and radiative cooling, and 
that it reaches a saturated value like seen in prior MAD simulations. 



Figure 3. Time and 0-averaged quantities are plotted for MADfHR. Rest 
mass gas flow is traced by black streamlines. The colored (green, black, 
and blue) thick lines correspond to time-averages of quantities Y = 
{u'\ut,b^/po], respectively, at values [Y], = {0,-1,1|. While near the 
BH the flow has [fo^/polr ^ I as averaged directly, the dense inflow has 
[b^/po]t ^ 1 at all other radii. Streamlines nearest the polar axes are cal¬ 
culated using velocity structure five cells away from the axis for better vi¬ 
sualization of the qualitative geometry of the jet flow. The flow exhibits 
equatorial asymmetry over many inflow times as found in thicker MAD 
simulations. The red magnetic field line contours show how the magnetic 
flux threads the black hole and inner disk. 

Fig. 0 shows the time and 0-averaged rest-mass velocity 
streamlines as well as several key boundaries associated with the 
accretion flow. This includes showing the jet boundary and the 
flow stagnation surface (where material has zero radial velocity, 
moves outward at higher latitudes, and moves inward at lower lat¬ 
itudes). Most importantly, the plot indicates the time-0-averaged 
structure of the magnetic field once the inner disk and BH have 
reached magnetic flux saturation. As with the thicker MAD model 
plotted similarly in Fig. 6 of [McKinney et aL]j2012[ >, the time-(/i- 
averaged steady-state field structure, as well as the evolution we 
see in movies, demonstrate ordered poloidal magnetic flux is easily 
transported through the disk from large radii. 


3.2 Long-term evolution of MAD state 

In order to study the thin MAD out to large radii to ensure an 
accurate quasi-steady radiative efficiency, we next considered the 
disk having an initially nearly MAD-level of magnetic field. Fig. 
[T] and Fig. show the initial and evolved quasi-equilibrium state 
of our high-resolution model that is initially nearly MAD inside 
r = 30rg. As the accretion flow evolves from its initial state, rota¬ 
tion drives outward angular momentum transport by turbulent mag¬ 
netic stresses via the MRI as well as by a wind. Mass, energy, an¬ 
gular momentum, and magnetic flux are accreted and ejected, and 
by ~ 5000rj/c magnetic stresses balance incoming gas forces so 
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Figure 4. Evolved snapshot of the MADiHR model at f » 70, OOOr^/c, with identical layout to Fig.j^ The MAD state is highly non-axisymmetric and clumpy 
with a weak (t/j ~ 1 %) jet. The rest-mass accretion rate and magnetic flux on the horizon reach the MAD state quickly, which was as desired so that we can 
study a quasi-steady MAD out to large radii over long time periods. The dimensionless magnetic flux on the horizon is in agreement with our MADfHR model, 
showing a disk initial condition with a relatively stronger poloidal field strength still reaches the same MAD state near the BH. 


the accretion is MAD out to r ~ 30rj. By the end of the simulation, 
t ~ 70000rj/c for MADiHR, the MAD extends to r ~ 35rj. 

The evolved state is evaluated using the period t = 
[30000,70000]r^/c for time-averaging quantities. A slightly longer 
period t = [10000,70000]rj/c is used for the time-averaged func¬ 
tions of radius in Fig.j^to minimize the effects of long-term trends. 
The time-0-averaged structure is show in Fig.|^and shows qualita¬ 
tively similar behavior to MADfHR. 

The gas-t-radiative efficiency rj ~ 20%, and specific angular 
momenfum j » -0.4 are constant within the MAD/non-MAD tran¬ 
sition r ~ 35rg, indicating a quasi-steady MAD state has been 
reached out to a relatively large radius. 

The disk thickness at r ~ 20r^ is regulated by cooling to be 
H/R X 0.11 for n=l in Eq. H/R « 0.14 for n - 2, density- 
weighted normalized H/R x 0.13 ( |Penna et al.|2012| >, and thermal 
thickness of h = arctanfcs/frot) ~ 0.13. An accurate measurement 
of H/R is complicated by the fact that magnetic compression leads 
to a geometric H/R x 0.05 near the BH, mass is launched into 
a broad wind, and vertical disk fluctuations are the same order as 
H/R. For the 0-averaged density profile and including only bound 
material (i.e. —u,(p + eg^s + P% + '2-Pb)/P < 1), valid when comparing 
with NT, we get H/R x 0.10, which is our target value we desired 
for the disk proper. Fig. shows the time evolution of H/R at four 
radii for n = 1. 

Our thin MAD achieves a magnetic flux on the BH of Th ~ 5. 


The stability of Th with time and the presence of plentiful magnetic 
flux remaining in the disk (see Fig. lower panel) indicate that 
the disk has achieved a balance between magnetic forces and gas 
forces with a strong magnetic field (T » 1). Thus one expects 
maximum magnetically-induced deviations from the standard NT 
thin disk solution. 

Large temporal deviations are driven by magnetic flux repeat¬ 
edly building up on the BH horizon to a level the disk cannot sus¬ 
tain. As soon as a weak point in the density exists caused by a non- 
axisymmetric mode, then there is a sudden re-emergence of mag¬ 
netic flux back into the disk, a process that creates large low-density 
bubbles in the form of a magnetic prominence. The magnetic flux 
that is ejected is then processed by turbulence by the outer parts of 
the disk until being accreted again at which point the process re¬ 
peats. Fig.j^shows one of the most extreme magnetic prominences 
occurring in the MADiHR simulation. 

3.3 Radiative Efficiency 

We find that nearly as much radiation is released near and inside the 
ISCO as in the entire rest of the disk, so we calculate the radiative 
efficiency assuming 100% photon capture within a choice among 
three radii: the horizon, the photon orbit radius r x 2.4, and the 
critical impact parameter radius, r » 4.1, for which geodesics with 
any smaller impact parameter intersect the horizon. An estimate of 
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Figure 5. Time-averaged quantities as a function of radius. From top to bot¬ 
tom, panels ai‘e: Total mass accretion rate (M, blue line), energy efficiency 
(;;), specific angular momentum ( 7 ), magnetic flux (T), and disk thickness 
(////?, n=l in Eq.[^ blue line). Plots for t\ and j show total gas-i-radiative 
(black solid), electromagnetic (green), matter (red), and radiation (black 
dashed) components. Remaining plots show their initial value (magenta dot- 
dash). The plot of T shows that for the total BH+disk contribution (blue). 
The constancy of M, rj, j show where the flow has achieved inflow equilib¬ 
rium. T shows the MAD to non-MAD transition at r ~ 35rg, while through¬ 
out the disk H/R w 0.1 (bound gas), as maintained by cooling. The primary 
result of this plot is the efficiency panel showing the radiative contribution 
of rji- ~ 15% with most of emission beyond the standard NT disk model 
coming from quite close to the photon orbit. 


photon capture using ray-tracing from a thin disk (see |Beckwith| 
|et al.|20Q6| t gives as much as 30% escape fraction from the pho¬ 
ton orbit for equivalent BH spin, and more than 70% for r » 4.1. 
This suggests using the photon orbit may give the best, simplest 
approximation. 

The radiative efficiency is t/i- -- 15% for MADiHR (< 0.5% 
lower for MADiMR) when including only radiation from outside 
the photon orbit. This is about twice larger than the rj,. ^ 8% effi¬ 
ciency in standard NT thin disk theory for BH spin ajM ^ 0.5. The 
impact parameter estimate is an extreme, albeit unphysical, limit 
giving efficiency rj^ ~ 9.4%. Also unphysical, including all radi¬ 
ation released outside the horizon gives rf^ ~ 19%. Unbound gas 
only contributes ~ 1% to t/j-, a smaller relative contribution com¬ 
pared to prior simulations ^Penna et al.|2010[|Zhu et al.|2012^ . 



Figure 6. Shows temporal evolution of a few quantities for our MADiHR 
model. The disk scale height HjR is shown in the top panel. We show 
the temporally-smoothed HjR as plotted for radii Risco (solid-black), 15rg 
(short-dashed), 30r^ (dots), and (long-dashed). The red line is the un¬ 
smoothed evolution at 15rg to demonstrate the level of variability. The HjR 
value is regulated by cooling to be HjR ~ 0.1, which is what is achieved 
for radii just outside the ISCO and out to the inflow equilibrium radius of 
r ~ 35rg. The middle panel shows the viscosity parameter spatially- 
averaged over the disk at a radius lOr^. The relatively high ~ 0.5 is 
indicative of MADs, implying efficient transport of mass and angular mo¬ 
mentum. The third panel shows the ratio of the total magnetic flux threading 
the BH to the total magnetic flux available at the start of the simulation in¬ 
side the radius out to which the simulation ultimately reaches at least one 
inflow time. This shows that the horizon contains much more magnetic flux 
than is available, so the quasi-steady state has a horizon-threading flux that 
has reached saturation. 


Thus, our thin MAD with HjR ^ 0.05-0.1 and a radiative 
efficiency of //r ~ 15%, demonstrates an 80% deviation from NT. 
Prior GRMHD simulations of thin disks have had HjR 0.1 with 
6 % deviations jNoble et al.|2009'l t, HjR ^ 0.06 with 6% to 10% 
deviations jNoble et al.|2011) , and HjR ^ 0.07-0.13 with 4% to 
5% deviations ( |Penna et al.|2010^ . So despite having similar HjR, 
our thin MAD has vastly higher deviations from NT. 


3.4 Jet Efficiency 

Some of the horizon magnetic flux reaches into the polar regions 
leading to a jet efficiency of zyj « 1% at r = 50r^ for both MADiHR 
and MADfHR. A caveat is that the jet within E jpo > 1 at r = 
is only resolved by about 7 grid cells. So we check the jet effi¬ 
ciency and resolution on the horizon, where the jet is resolved with 
more than 12 grid points per hemisphere and still has a jet power 
of 77 ~ 1%. These horizon jet grid cells are sufficient resolution to 
capture the force-free behavior of the jet near the polar axes, ac¬ 
cording to prior resolution studies ( [McKinney & Gammie|2004| ). 
So the jet efficiency measurement is robust to resolution changes, 
slight model changes, and what radius it is measured at. 

The BH spin and disk rotation drive a magnetized wind with 
efficiency == 4%, though at least an additional 1% wind effi¬ 
ciency may be lost to radiation via the very efficient cooling we use 
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Figure 7. Time and 0-averaged quantities are plotted for MADiHR. Iden¬ 
tical layout as Fig. The disk is magnetically arrested out to a transition 
radius of ~ 35rg, which leads to a short inflow time-scale compared to the 
MADfHR disk at those radii, so the time-averaged velocity fields and other 
disk structure is smoother and more symmetric. 


as well as lost to the comparatively low resolution in the jet and 
highly magnetized wind regions. 

The jet and wind contribute only a bit to the total rj x 20%. 
However, the jet and wind carry nearly as much angular momentum 
outward as the dense Keplerian disk carries inward. This results in 
the BH spin-up parameter s = -j - 2(a/M)(l — rj) x -0.4, i.e. 
somewhat spinning down ( [McKinney et al.|2012^ . 

Further, the jet efficiency as a function of HjR is constrained 
by our thin MAD result. Magnetic compression leads to varying 
geometric HjR vs. radius, so the more constant th ermal thickness 
li = arctan (cjvmt) is used i McKinney et al.|2012 i. Our thin MAD 
models have h x 0.13 and older MADs had h x 0.5 (e.g. with 
geometric HjR x 0.3 at r = 20?-j) orhx 1.5. 

Prior thick MAD simulations found that the saturated mag¬ 
netic flux on the BH given by T is roughly constant vs. BH spin. 
That data showed a slight tendency to peak at ajM ~ 0.5, but this 
is probably intrinsic to the initial disk conditions or due to intrinsic 
errors that suggest ajM ~ 0 could be the peak as well [McKinn^ 
jet al.|201^|Tchekhovskoy et al.|2012l (. The thick MAD simulations 
of moderate thickness show a clear asymmetry in the T vs. spin for 
positive vs. negative spin that is consistent with the asymmetry in 
jet efficiency vs. spin. So we do not try to fit the peak in T being at 
ajM ~ 0.5, but we capture the asymmetry. We also capture the fact 
that this spin asymmetry is seen to go away as the disc becomes 
even thicker. 

Combining the different thick MAD models with our new thin 
MAD model, a rough fit of the jet power is given by 


. , 0.3m„ 

■ 400%mf, 1 -H-^ 

" ’ I+ 2k^ 




(13) 


where wh = aim , which assumes steady activation of a|Blandford| 


& Znajek (1977 i jet. The rjj oc scaling is required (compared to 
oc h) to reach our thin MAD with rjj x 1%. That r/j oc {HjR)^ is 


consistent with S7.1 ofjPenna et al.|([20T0l( who analytically found 


T oc (HIR), and with 


Meier 


120011 who argued that oc {HIRy. 


This formula is roughly consistent with what is expected based 
upon the BZ formula r/j oc Models with h x 0.5 and ajM = 

0.5 gave t/j ~ 20% and T x: 13. Our thin MAD has T 5, so the 
BZ formula suggests we should get rjj x 3% while we got rjj x 1 %. 
Why the factor of three difference? 


First, Y measured on the BH is not relevant to the BZ jet. One 
cannot use Yh with r/j oc Y^, because in general (and as found in 
thicker MAD models) magnetic flux in the jet is a fraction of Yh 
and the jet only partially covers the rotating BH’s surface. Instead, 
the magnetic flux passing through the jet is what is relevant to the 
jet power. The h x 0.5 model had Yj 10, while our model has 
Yj » 4. These are significantly different from the magnetic flux 
over the entire horizon. This still implies we should have rjj x 3% 
for consistency with the BZ expectation. 

Second, the BZ formula fails to account for the disk covering 
a significant part of the potential power output [McKinney|2005| (. 
The disk at the equator gets in the way, even for MADs (especially 
at lower spin when less magnetic compression occurs). However, 
most of the BZ power would be at the equator where the Poynting 
flux peaks. The disk chops-off this peaked part of the Poynting flux 
for the jet, giving some of it to the wind. Some of that BZ power 
is never produced in the first place due to the heavy mass-loading 
by the disc. A better model of that disk part that connects to the 
horizon is the [Gammiej j 1 999) model of a thin magnetized disc. In 
that case, Y locally in the disk is much less than over the entire 
horizon, so the |Gammie| ( |1999t model predicts much less than BZ 
for extraction of energy. One cannot just fit t/j vs. Yh and expect it 
to match BZ. 


Third, we have defined the jet arbitrarily as where K jpo > 1. 
The wind is defined as the rest of the non-radiative flux that adds 
up to the total efficiency. For our thin MAD, this is 4% (i.e. 20%- 
1%-15%). The wind efficiency is partially directly due to the disk 
rotation, but some of the wind efficiency derives from magnetic 
flux threading the BH that passes into the disk leading to a BH- 
spin-driven wind. Some of this BH-spin-driven wind absorbs what 
would have otherwise (without a disk) become the jet and can ac¬ 
count for another missing 2%, giving a total BZ-like efficiency 
matching expectations. If this were the sole effect in comparison 
to the thicker MAD models, this would suggest that thin MADs 
may have more powerful winds per unit jet power relative to thick 
MADs. 


Fourth, the jet efficiency for the h x 0.5 ajM = 0.5 model 
is a peak among all the otherwise identical models with different 
BH spins. Given how much the results vary for slight changes in 
the model setup (like initial magnetic field strength), the jyj 20% 
could easily have been t/j x: 10% within systematic errors. Then 
there would be no discrepancy between our formula and the BZ 
approximation. 

So our Eq. ( |13[ is not supposed to be derived from BZ because 
the BZ-driven process creates both a “jef’ and a “wind”, while we 
have defined the jet by that part which would be capable of becom¬ 
ing relativistic. So in that sense, our formula is useful in application 
to real jets instead of a test of whether the theoretical BZ extrac¬ 
tion efficiency applies to the simulation data. However, only direct 
observational comparisons without ad hoc definitions will test the 
validity of these simulations being applicable to real systems. The 
fitting formula is just a useful guide. 
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Figure 8. Plotted as in the top frames of Fig.j^is a snapshot of MADiHR at t = 31568rj/c when the disk experiences an extreme eruption of magnetic flux 
from the BH horizon. A large magnetic flux tube creates a low-density cavity in the disk before it is parfially blended with disk material by the magnetic 
Rayleigh-Taylor instability. Such prominences occur frequently in the MAD state. 


3.5 Resolution 


We test the resolution of our simulations using convergence quality 
factors for the MRI (Qs^.mriI number of grid cells per MRI wave¬ 
length in each direction) and turbulence (Qre(p,corr', number of cells 


per turbulence correlation length in each direction) ( 

Shiokawa et al. 

2012 

[Hawley et al. 2011 

IMcKinney et al.||2012 

Hawley et al. 

2013 

[McKinney et al.|20 

3 2014b. For our high resolution runs. 


MADiHR (MADfHR), at r = 12, 30rg centered on the equatorial 
flow, we have Qa^mri = 160,55 (Qb.mri = 88,18) with a weighting 
of and we have = 120,100 (2 ,^,mri = 79,77) with 


a weight of p. Weighted averages are computed for the numerator 
and denominator separately before taking the ratio to get Q. The 
purely density-weighted values represent a stricter constraint on 
whether we have sufficient resolution (because MADs have strong 
fields mostly outside the dense region, and the external stresses are 
important to whether transport occurs), and these are plotted as a 
function of radius in Fig.|^ Resolution of rest-mass and magnetic 
energy densities are evaluated at four radii r = {r/;/rj,4,8,30)rj. 
For MADiHR, rest-mass correlation values are, at each radius re¬ 
spectively, Qr,corr ~ 6,10,14,14, Qg^corr ~ 40,35,31,34, and 
Q,p,corr ~ 6,7, 8,11. Correlation values for magnetic energy den¬ 
sity are, similarly, ~ 6,12,17,15, Qg^corr ~ 34,33,29,36, 

and Q^,corr ~ 8,9,12,16. 

For MADfHR, calculated for the disk once the BH has reached 
the MAD state, rest-mass correlation values are, at each radius 
respectively, Qr,corr ~ 6,10,14,14, Qg^corr ~ 3,29,28,32, and 
Q,p,corr ~ 7,9,11,13. Correlation values for magnetic energy den¬ 
sity are, similarly, Qr,corr ~ 6,12,16,15, Qg^corr ~ 5,28,30,31, and 

Q^.corr~ 1,9,12, \A.’ 

These MRI and turbulent quality factors are high enough to 
ensure the reported simulations are resolving all the necessary tur¬ 
bulent physics and are likely converged. We also directly conver¬ 
gence tested our results by performing simulations of each setup 
with three resolutions. The qualitative structure of the disk in the 
MR and HR simulations are the same, but the lower resolution in 


MR allowed us to run these simulations for significantly longer. 
The nature of magnetic flux accumulation is the same, but the 
longer MR simulations become magnetically arrested to larger ra¬ 
dius, r ~ 18r^ in MADfMR compared to only r ~ 15rj in 
MADfHR. A difference of only ~ 0.5% in radiative efficiency be¬ 
tween MADiMR and MADiHR indicates very good convergence 
in our quantity of interest. We find only small (typically ~ 10% 
or less) differences between our MR and HR simulations for the 
other quantities reported which characterize the accretion flow. For 
instance, the jet power differs by only ~ 50% befween MADiMR 
and MADiHR, unsurprising since T 5 for both. Low resolution 
‘testing’ simulations were largely in qualitative agreement, but with 
some 2’s indicating lack of convergence. 

Figs.|^and|^show values of the viscosity parameter, dom¬ 
inated by the magnetic stress. We find Ub ^ 0.5, as consistent 
with resolved disks in prior MAD simulations. The effective vis¬ 
cosity Q!b,eff. corresponding to a measurement of the actual inflow 
rale indicated by an a parameter, is much larger than at, due to 
enhanced inflow caused by external torques driven by low-density 
highly-magnetized material that introduces large-scale stresses on 
the more dense inflow. 


3.6 Power-law Fits for Radial Dependence 

We now consider radial power-law fits (of the form / = /ofr/ro)") 
for quantities associated with the disk structure and flow, plotted in 
Fi g.|10| for MADiH R. The fits to disk quantities are performed as 
in [McKinney et aL|^2012^ between radii r = 12rg and r = 30rg 
for both our high resolution models. The profiles for MADiHR 
(MADfHR) arepo oc ^-0-4±o.i ^-o. 8 ±o. 2 )^ ^-i.5±o.i 

\Vr\ oc (oc OC ^-"38*0.01 ^-0.49±0.01) ^ 

\b,.\ oc 1^1 ^ ^-1.240.05 ^-0.89±0.02) 

The MADiHR fits are more representative of a steady-state 
MAD solution, because the disk in MADiHR is actually magnet¬ 
ically arrested over all radii used for fitting and the model is in 
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Figure 9. MRI quality factors indicating time-(^ averaged resolution, vis¬ 
cosity of the disk, and vertical MRI half-wavelengths per disk scale height 
(Sd, blue line) are plotted as a function of radius. Q factors indicated the 
number of grid cells per critical MRI wavelength, averaged with weight 
(solid lines) and p (short dashed) across the disk where jp < 0.5. 
The third panel shows the viscosity parameter o-b (solid line) and the ef¬ 
fective viscosity calculated via the mass accretion rate (short dashed). The 
magnetic evolution of the disk is well resolved and leads to a high viscosity 
as measured by local stresses, a^, and the rate of mass inflow O'b.eff ■ The ma¬ 
genta dot-dash line represents the value of 5 d at t=0 and the time-average 
is blue-solid. 



Figure 10. Several quantities characterizing the accretion flow in MADiHR 
are plotted as a function of radius. The top panel shows the density profile po 
(solid), the internal gas energy Ug (dashed), and magnetic energy density Ub 
(dotted). The inner disk, that has reached the MAD state, has Ug ~ Ub indi¬ 
cating equipartition is achieved even within the dense parts of the disk. The 
middle panel includes the disk radial velocity, azimuthal velocity, and Ke- 
plerian velocity (Vi-, solid; black-dashed; vk, blue-dashed, respectively). 
The final panel shows the radial profiles of each component of the comov¬ 
ing magnetic field, in solid, dashed, and dotted lines respectively. The 
density is shallower than thicker disk solutions, the disk is fairly Keplerian, 
and the magnetic field falls off as expected for accretion disks. 


inflow equilibrium out to large radii. Our thin disk has a somewhat 
shallower radial density profile than prior MAD models, but are 
otherwise consistent. To obtain a density profile accurate to large 
radii, much longer-run simulations would be required. The overall 
fits are consistent with pg oc and \b\ oc expected for accre¬ 
tion disks, while the shallow density profile could be connected 
with the stronger wind (relative to the jet) compared to thicker 
MAD simulations. 


4 SUMMARY AND CONCLUSIONS 


We have performed fully 3D GRMHD simulations of radiatively 
efficient thin accretion disks seeking to maximize the magnetically- 
induced deviations from NT thin disk theory. We first demonstrated 
that even in thin disks, after a time 108,000rg/c and out to r ~ 
17rg, magnetic flux readily adverts inward and builds-up to a MAD 
level (for which accumulation results in magnetic forces pushing 
out balancing against gas forces pushing in). This occurs despite 
possible outward magnetic diffusion through the disk. 

Our simulation that is MAD out to r ~ 35rg after 70, OOOrg/c 
has a radiative efficiency of x: 15%. This is 80% higher than 
the standard NT thin disk value of pnt ~ 8% and is as if the disk 
were a standard NT thin disk but with ajM x 0.9. BH X-ray bina¬ 
ries have LIL^ii ~ 0.1 in the high-soft state ( [Kulkami et al.|20TT| ), 
for which the NT model gives HjR < 0.02 near the BH. Prior non- 
MAD simulations suggest deviations scale with HjR jPenna et al.| 
[ 20 T 0 I 1 , so typical thin MAD deviations might be ~ 20%. Such a 
magnitude of deviations could slightly impact BH spin measure¬ 
ments, such as in [Noble et al.| ( [2011^ who found 6% deviations 
lead to a BH with a/M = 0 giving a spectrum like NT but with 
ajM X 0.3. |Zhu et al.|j^l2^ found that emission from near the 
ISCO is mostly high-energy, but more physics (e.g., dissipation- 
radiation energy balance, optical depth effects, BH photon capture) 
is required for BH spin fitting ( [Kulkami et al.|2011^ . We plan to 
obtain spectra via post-processing jZhu et al.|2015^ . Also, clumpy 
rest-mass distribution in thin MAD flows may help to explain spec¬ 
tra from AGN 1 Dexter & Agol|20lT) and BH X-ray binaries jDex-[ 
[ter & Quataert|2012p . 

It is not a new idea that X-ray binaries in the low-hard (LH) 
state supply enough magnetic flux to the BH to enable BZ-type 
powering of an observable jet by being geometrically thick, and 
that the jet quenches when the disk collapses to become thin in 
the high-soft state. For instance, plunging material inside the ISCO 
may enhance the flux on the BH in the low-hard state only, even 
in sub-MAD disks ^Reynolds et al.|[2006^ . However, these ideas 
generally propose transition mechanisms in which the jet is either 
powered or not, depending on disk state. Our model does not re¬ 
quire a special mechanism for the LH state. 

Radiatively inefficient flows have H/R x 0.4 (jNarayan et al.| 
[2012[ (. Our fitting formula Eq. suggests that jet power in 
the MAD state can drop by > (0.4/0.02)^ ~ 400 (or drop by 
> (0.4/0.01)2 _ igQQ j-Qj. 4 u 1957+11 if 0 .O 6 ) when 

undergoing hard-to-soft state transitions. This implies that the ob¬ 
served jet quenching in BH X-ray binaries jRussell et al.|201 1) can 
occur despite the presence of a large-scale MAD-level magnetic 
field in both disk states. Such hard-to-soft transitions may also be 
applicable to tidal disruption events at late times ([Tchekhovskoy[ 
[et al.|2014^ . The MAD state could be necessary to explain power¬ 
ful jet systems ( [Zamaninasab et al.|2014[ (, while the thermodynam¬ 
ically determined disk thickness may also play an important role 
in setting just how much magnetic flux a MAD has. Even thinner 
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disks than studied here would be implied to have T oc HjR and 
so could have much weaker magnetic flux and much weaker jets 
rather than being independent of the disk state. 
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